#include <fstream>
#include <iostream>

float wsp_sprez(bool napr,int f,float e,float v)
{
 if(napr)
 {
  switch (f)
  {
   case 1: return (e)/(1-(v*v)); break;
   case 2: return (e*v)/(1-(v*v)); break;
   case 3: return (e)/(2*(1+v)); break;
   default: cout << "Cos zrypales przy sprzezystosci naprezenia" << endl; return 1; // THROWNAC
  }
 }

 if(!napr)
 {
  switch (f)
  {
   case 1: return (e*(1-v))/((1+v)*(1-(2*v))); break;
   case 2: return (e*v)/((1+v)*(1-(2*v))); break;
   case 3: return (e)/(2*(1+v)); break;
   default: cout << "Cos zrypales przy sprzezystosci odksztalcenia" << endl; return 1; // THROWNAC
  }
 }
 return 1;
}

float LiczPole(punkt a,punkt b, punkt c)
{
  return 0.5*((a.x*(b.y-c.y))+(b.x*(c.y-a.y))+(c.x*(a.y-b.y)));
}

float C2(int i, punkt a, punkt b, punkt c)
{
  float A=LiczPole(a,b,c);
  if(i==1)return ((b.y-c.y)/(2*A));
  if(i==2)return ((c.y-a.y)/(2*A));
  if(i==3)return ((a.y-b.y)/(2*A));
  exit(-1);  //tu tez trzeba cos thrownac
}

float C3(int i, punkt a, punkt b, punkt c)
{
 float A=LiczPole(a,b,c);
 if(i==1)return ((c.x-b.x)/(2*A));
 if(i==2)return ((a.x-c.x)/(2*A));
 if(i==3)return ((b.x-a.x)/(2*A));
 exit(-1);  //tu tez trzeba cos thrownac
}
vector<vector<float> > KijE(siatka & s,stale per,int i,int j)
{
 vector<float> v(2);
 vector<vector<float> > mac(2,v);
/* float A=LiczPole();
 mac[0][0]=(wsp_sprez(true,1,per.E,per.V)*C2(i,)*C2(j,)+wsp_sprez(true,3,)*C3(i,)*C3(j,));
 mac[0][1]=(wsp_sprez(true,2,per.E,per.V)*C2(i,)*C3(j,)+wsp_sprez(true,3,)*C3(i,)*C2(j,));
 mac[1][0]=(wsp_sprez(true,2,per.E,per.V)*C3(i,)*C2(j,)+wsp_sprez(true,3,)*C2(i,)*C3(j,));
 mac[1][1]=(wsp_sprez(true,1,per.E,per.V)*C3(i,)*C3(j,)+wsp_sprez(true,3,)*C2(i,)*C2(j,));
*/
 return mac; 
}

//Wspol C
/*float * siatka::wspol(int nrele)
{
 float A;//Pole pojedynczego elementu
 float * tab=new float [7];
 float x1,x2,x3,y1,y2,y3;
 int tempinw;//chwilowy nr wezla tworzacego
 tempinw=kon[nrele].wt1;
 x1=wsp[tempinw-1].x;
 y1=wsp[tempinw-1].y;
 tempinw=kon[nrele].wt2;
 x2=wsp[tempinw-1].x;
 y2=wsp[tempinw-1].y;
 tempinw=kon[nrele].wt3;
 x3=wsp[tempinw-1].x;
 y3=wsp[tempinw-1].y;
 A=0.5*(x1*(y2-y3)+(x2*(y3-y1))+(x3*(y1-y2)));
 //C21
 tab[0]=(y2-y3)/(2*A);
 //C22
 tab[1]=(y3-y1)/(2*A);
 //C23
 tab[2]=(y1-y2)/(2*A);
 //C31
 tab[3]=(x3-x2)/(2*A);
 //C32
 tab[4]=(x1-x3)/(2*A);
 //C33
 tab[5]=(x2-x1)/(2*A);
 //Pole elementu do uzycie w innych wzroach
 tab[6]=A;
return tab; 
}*/
/*
vector<vector<float> > MacierzGlobalna(siatka & s)
{
 
}

//zwraca macierz
float ** MacGlob(siatka & s)
{
 int roz=s.liczba_wezlow*2;
 float ** temp=new float* [roz];
 for(int i=0;i<roz;i++)
  {
   temp[i]=new float [roz];
  }
 for(int i=0;i<roz;i++)
  {
   for(int j=0;j<roz;j++)
    {
     temp[i][j]=0;
    }
  }
 for(int i=0;i<roz;i++)
  {
   int nrelem=s.koneksje[i].nr;
   float *tabWSP=wspol(nrelem);
   int tabwez[6];
   //i11
   tabwez[0]=(kon[nrelem].wt1*2)-2;
   //i12
   tabwez[1]=(kon[nrelem].wt1*2)-1;
   //i21
   tabwez[2]=(kon[nrelem].wt2*2)-2;
   //i22
   tabwez[3]=(kon[nrelem].wt2*2)-1;
   //i31
   tabwez[4]=(kon[nrelem].wt3*2)-2; 
   //i32
   tabwez[5]=(kon[nrelem].wt3*2)-1;
  int y=0;
  for(int x=0;x<3;x++)
  {
   int w=0;
   for(int z=0;z<3;z++)
    {
     temp[tabwez[y]][tabwez[w]]+=LiczPole()*((f1*tabWSP[x]*tabWSP[z])+(f3*tabWSP[3+x]*tabWSP[3+z]));
     temp[tabwez[y]][tabwez[w+1]]+=LiczPole()*((f2*tabWSP[x]*tabWSP[3+z])+(f3*tabWSP[3+x]*tabWSP[z]));
     temp[tabwez[y+1]][tabwez[w]]+=LiczPole()*((f2*tabWSP[3+x]*tabWSP[z])+(f3*tabWSP[x]*tabWSP[3+z]));
     temp[tabwez[y+1]][tabwez[w+1]]+=LiczPole()*((f1*tabWSP[3+x]*tabWSP[3+z])+(f3*tabWSP[x]*tabWSP[z]));
    w+=2;
    }
  y+=2;
  }  
  delete [] tabWSP;
  }
  ofstream wyj;
  wyj.precision(6);
  wyj.flags(ios::scientific);
  wyj.open("test2");
 for(int i=0;i<(roz);i++)
  {
   for(int j=0;j<(roz);j++)
    {
     wyj<<temp[i][j]<<"\t";
    }
   wyj<<"\n";
  }
 wyj.close(); 
return temp; 
}
*/
